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ABSTRACT 

We report on Nuclear Spectroscopic Telescope Array (NuSTAR) hard X-ray observations of the 
young rotation-powered radio pulsar PSR B1509—59 in the supernova remnant MSH 15—5 2. We 
confirm the previously reported curvature in the hard X-ray spectrum, showing that a log parabolic 
model provides a statistically superior fit to the spectrum compared with the standard power law. The 
log parabolic model describes the NuSTAR data, as well as previously published y-ray data obtained 
with COMPTEL and AGILE, all together spanning 3 keV through 500 MeV. Our spectral modelling 
allows us to constrain the peak of the broadband high energy spectrum to be at 2.6±0.8 MeV, an 
improvement of nearly an order of magnitude in precision over previous measurements. In addition, 
we calculate NuSTAR. spectra in 26 pulse phase bins and confirm previously reported variations of 
photon indices with phase. Finally, we measure the pulsed fraction of PSR B1509—58 in the hard 
X-ray energy band for the first time. Using the energy resolved pulsed fraction results, we estimate 
that the pulsar’s off-pulse emission has a photon index value between 1.26 and 1.96. Our results 
support a model in which the pulsar’s lack of GeV emission is due to viewing geometry, with the 
X-rays originating from synchrotron emission from secondary pairs in the magnetosphere. 

Subject headings: pulsars: individual (PSR B1509—58) - stars: neutron - X-rays: stars 


1. INTRODUCTION 

PSR B1509—58 is a young, energetic rotation-powered 
pulsar discovered i n the soft X-r ay ba nd usin g the Ein¬ 
stein Observatory (Seward & Harnden 1982[) and soon 
afterw ards detected in the radio band (Manchester et al. 
1982). It has period ~150 ms, period derivative ~1.5 
xl0 _12 s s _1 and a high spin-down luminosity of E = 
1.7 x 10 37 crg s _1 . Assuming the standard magnetic 
dipole model, its characteristic spin-down age is ~1600 
years and inferred surface magnetic field at the equator 
is ~ 1.5 x 10 13 G, higher than the lowest e stimat ion of 
£?dip — 6 x 10 12 G for a magnetar (Rea et al.||2013). The 
pulsar lies in the centre of the supernova remnant MSH 
15—5,9, which has a complex structure with thermal and 
non-the rmal emissiom whic h is largely powered by the 
pulsar (Tamura et al. 1996). Recently, NuSTAR obser¬ 
vations ol the M SH 15—59 and its pulsar wind nebula 
(PWN) region by An et al. (20141 revealed clear evidence 
for synchrotron burnoff m the PWN, spectral softening 
with distance from the pulsar, and an interesting shell¬ 
like structure in the Ah map. 

The aforementioned NuSTAR observations also offer a 
new window into the hard X-ray emission of the central 
pulsar itself, the subject of this work. PSR B1509—58 
is one of the brightest rotation-powered X-ray pulsars 
in the sky and has been observed by ma ny X-ray and 
soft y -ray tele scopes, for example , Ginga (Kawai et al. 
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dKuiper et al.||1999|), AGILE (j Pilia et al.||2010|), Bep¬ 
poSAX (|Cusumano et al.||2001 1, and Fermi LAI’ (den 
Hartog et al. in preparation). 'Thanks to NuSTARA wide 
energy range in the hard X-ray band and unique hard 
X-ray focusing ability therein, NuSTAR observations of 
PSR B1509—58 can shed new light on the pulsar’s hard 
X-ray spectrum. 

The pulsar’s spectrum in the hard X-ray band has long 
been of interest, in part because the source’s large flux 
enables detailed studies not possible in most rotation- 
powered puls ars. Based on BeppoSA X and COMPTEL 
observations, Cusumano et al. (20011 reported that the 
spectral energy distribution of PSR B1509—58 is well 
represented by a logarithmic parabolic function char¬ 
acterized by a linear dependence of the spectral slope 
upon the logarithm of energy from ~1 keV to ~30 MeV, 
rather than the simple power-law model so ubiquitously 
assumed for rotation-powered pulsar spectra. The spec¬ 
tral bending in PSR B1509—58 was reported to be de¬ 
tectable also in the keV energy range. If correct, this 
has important implications for the high energy emission 
mechanism. However, BeppoSAX was an X-ray astron¬ 
omy satellite that consisted of three instruments oper- 
ating in different energy bands ( Boella et al.|[l9 97), and 
slight systematic cross-calibration discrepancies among 
these three instruments have been reported by |Kirsch 
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et al. (2005|. Hence confirmation of the claim of spectral 
curvature from BeppoSAX data is important to obtain, 
particularly with a mission like NuSTAR which covers 
the full 3-79 keV band with a single instrument. This, 
when combined with recent soft 7-ray observations, can 
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in principle provide vastly improved constraints on the 
full high-energy spectrum, which could help better deter¬ 
mine the peak spectral energy and ultimately address the 
curious low-energy cutoff s een in the high-en ergy 7 -ray 
spectrum from this source (Abdo et al.|2010 1 . 

In this paper, we present timing and spectral analyses 
of PSR B1509—58 data obtained by NuSTAR as well as 
the Chandra X-ray Observatory. In Section 2] we de¬ 
scribe our observations and data reduction. In Section 
[3j we present our analysis and results. In Section [4] we 
discuss the results and in Section [5] provide a summary. 

2. OBSERVATIONS 

NuSTAR is the firs t focusing X-ray tele scope that op¬ 
erates above 10 keV ( Harrison et al.|2013 ). It consists of 
two co-aligned hard X-ray grazing incidence telescopes, 
which focus onto two independent solid-state focal plan 
modules, and produce two spectra in each observation. 
The angular resolution of the observatory is 18" FWHM, 
with a half power diameter of 58". Both telescopes oper¬ 
ate in the energy band from 3 to 79 keV and each provides 
an energy resolution of 400 eV FWHM at 10 keV and 0.9 
keV at 60 keV. The temporal resolution is 2 /is, and the 
timing accuracy is 1-2 ms. 

Four NuSTAR observations were made of the region 
containing PSR B1509—58 as well as its PWN in June 
and Au gust of 2013, wi th a total exposure time of 178.9 
ks (see An et al. 2014). In every observation, the two 
telescopies each produced a spectrum for a total of eight 
spectra. We also analyzed an archival Chandra observa¬ 
tion made with the High Resolution Camera (HRC) in 
June 2005, with an exposure time of 44.9 ks. The HRC is 
sensitive within 0.08-10 keV. Since the spectral response 
of HRC is extremely limited, we used the Chandra data 
only for the timing analysis below. Table [l] summarizes 
the observations. 

The NuSTAR data were processed with HEASOFT 
v6.15.1, nupipeline vl.3.1, and the Calibration 
Database (CALDB) files from 20131223. Standard Level 
2 cleaned events files were generated for further analysis. 
The Chandra data were reprocessed with chandra.repro 
in Cl AO 4.5 in order to use the most recent calibration 
files. 


3. DATA ANALYSIS AND RESULTS 
3.1. Timing Analysis 

From each NuSTAR observation, pulsar events were 
extracted in the 3-79 keV energy band using a circular 
region of radius 30" centered on the pulsar. We then ap¬ 
plied a correction to convert photon arrival times to the 
equivalent time at the solar system barycenter (Barycen- 
tric Dynamical Time) using the DE200 ephemeris. We 
verified that the pulse periods at each epoch were consis¬ 
tent with those predicted by the ephemeris of Livingstone 
& KaspilpMTl). 

To measure the properties of the X-ray pulsations of 
PSR B1509—58, the background must be carefully sub¬ 
tracted. The background consists of two components: 
detector background and PWN. For the NuSTAR obser¬ 
vations, the detector background was extracted with an 
aperture of radius 45" in a pulsar-free and PWN-free re¬ 
gion on the same detector chip with the pulsar. We also 
tried different background regions and found that all re¬ 
sults in our analysis were independent of the exact choice 



Figure 1. Pulse profiles from NuSTAR observations in various en¬ 
ergy bands: from bottom to top are energy band of 3-79 keV, 3-4 
keV, 4-5 keV, 5-6 keV, 6-7 keV, 7-8 keV, 8-9 keV, 9-10 keV, 10-11 
keV, 11-12 keV, 12-14 keV, 14-16 keV, 16-19 keV, 19-25 keV and 
25—79 keV, respectively. Phase zero was selected approximately 
around the point where the pulsation starts. The histograms each 
contain 64 bins. Values of the y-axis don’t have physical meaning. 
No significant change in shape was found in different energy band- 
passes. We used phas e 0.7 -1 (vertical dotted line) as the off-pulse 
component in Section |3.2| 


of region. The PWN, however, is important to consider 
because of NuSTAR’ s broad point spread function (PSF), 
which will be described later in this section. 

We produced pulse profiles using the NuSTAR. data in 
14 energy bandpasses from 3-79 keV (Fig. l}. Each pulse 
profile contained 64 bins and the number ofevents in each 
bin was at least 20. To obtain the pulsed component from 
a pulse profile, we determined the smallest photon count 
value among the 64 bins, used it as the off-pulse level, 
and subtracted this off-pulse value from each of the 64 
bins. 

In order to search for pulse profile evolution with en¬ 
ergy, we compared the shapes of the pulse profiles in the 
14 energy bands. We scaled them by dividing each with 
the total photon counts in that pulse profile. We chose 
the pulse profile of 3-79 keV as a reference, subtracted 
each of the 14 normalized profiles from it, and tested y 2 . 
No significant variation in shape was detected, and the 
residuals appeared to be random. The y 2 test showed 
that the probability of the pulse profiles being the same 
in shape with the reference is above 40% for 11 energy 
bands. The most suspicious one had a probability of 
0.24%, still not sufficiently small to indicate a disagree¬ 
ment in shape, considering the number of trials. Further¬ 
more, shifts in phase could indicate that emission regions 
vary with energy. To check this, we performed a cross 
correlation of each pulse profile with the reference, but 
no significant phase shift was seen for any profile. 

The pulsed fraction can be defined in several ways. 
Here we chose the RMS pulsed fraction, as well as an 
intuitive definition based on the area under the pulse . 
The RMS pulsed fraction is defined as (Dib et ah||20 09) 


PFrms = 
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Table 1 

Summary of Observations 


Observatory 

Obs. ID 

Start Epoch 
(MJD) 

Start Date 
(UT) 

Stop Date 
(UT) 

Exposure 

(ks) 

NuSTAR 

40024004002 

56450.9 

2013Jun07 

2013Jun08 

43.4 

NuSTAR 

40024002001 

56451.8 

2013Jun08 

2013Jun09 

42.6 

NuSTAR 

40024003001 

56452.6 

2013Jun09 

2013Junl0 

44.3 

NuSTAR 

40024001002 

56519.6 

2013Augl5 

2013Augl6 

48.6 

Chandra (HRC) 

5515 

53534.0 

2005Junl3 

2005Junl3 

44.9 



Energy (keV) 


Figure 2. RMS and area pulsed fractions in multiple energy bands 
from NuSTAR and Chandra observations. On the one hand, the 
NuSTAR results increase significantly between 3 to ~19 keV and 
then approach a plateau within 19-79 keV: the RMS pulsed fraction 
of the highest energy band is about 3 cr larger than that of the 
lowest one. On the other hand, the Chandra results do not show 
any variation between 0.08-10 keV. Moreover, the RMS pulsed 
fraction of Chandra is consistent with NuSTAR's results above 19 
keV within 0.2a. 

Here, a k =jf J2iLi r i cos(27n^), b k =jf J^iLi r i sin(27r0i), 
and (J ak , <Jb k are the uncertainties in a k and b k respec¬ 
tively. <f>i is the phase, N is the total number of phase 
bins, and k refers to the number of Fourier harmonics. 

The area pulsed fraction is 


- r min) 

y N r 

Lji— i '* 


PF — 

J r area — 


Here, r.\ is the count number in the <th phase bin, N is 
the total number of phase bins, and r min is the minimal 
count number in all bins. 

Figure |2] displays the RMS and area pulsed fractions 
measurecffor pulse profiles in the 14 energy bands. Both 
the RMS and area pulsed fractions increase with energy 
between 3 and ~19 keV and then approach a plateau at 
19-79 keV. The RMS of the highest energy band is about 
3er larger than that of the lowest band. 

To estimate the pulsed fraction that is free from PWN 
contamination (PWN-free), pulsed fractions were also 
measured from the Chandra HRC observation. The 
RMF of the Chandra HRC is complicated, so the energy- 
resolved HRC pulsed fractions (Fig. [2]) should be re¬ 
garded as as having been measured just in “soft” and 
“hard” bands, not the exact energy ranges. This does 
not change any result though, since we found the Chan¬ 
dra pulsed fraction value to be constant in these two 


Table 2 

Cross Normalization Factors of Each Spectrum 


Obs. ID 

Cross Normalization Factor 

40024004002A 

40024004002B 

40024002001A 

40024002001B 

40024003001A 

40024003001B 

40024001002A 

40024001002B 

i 

1.00 ±0.014 
0.789±g;°i| 

0.794 ± 0.013 

1.02 ±0.015 

1.04 ±0.015 

1.02 ±0.016 

1.03 ±0.016 


bands. For the Chandra observation, we can ignore back¬ 
ground contamination because of Chandra's narrow PSF. 
Specifically, in the Chandra observation within 0.08-10 
keV, there were 12100 events collected in the pulsar re¬ 
gion (i?=1.5") while the number of background events 
collected within an annular region with R in =2" and 
Rout= 4" was 374. Hence the background contamination 
in the Chandra pulsar region would be 70/12100, which is 
ignorable, and the PWN-free pulsed fraction value could 
be measured directly. 

We estimated the PWN contamination in 3-10 keV 
in the NuSTAR pulsar region using the Chandra HRC 
observation. In 0.08-10 keV, the Chandra observation 
had 8080 counts in an annular region with Ri n = 1.5 ,/ and 
Roiit=30". We used PIMM^/]to convert this into the 3- 
10 keV NuSTAR counts ass uming a p ower-law spectrum 
with a photon index of 1.7 (An et al.||2014|), and the re¬ 
sult was 12694 counts. In 3-10 keV, the NuSTAR. data 
had 37927 counts in the pulsar region, so the PWN con¬ 
tamination is roughly 12694/37927, not ignorable. This 
implied a pulsed fraction of ^66%, consistent within 3<r 
with our NuSTAR measurements in 3-4 keV, 4-5 keV, 
5-6 keV, 6-7 keV, 7-8 keV, 8-9 keV and 9-10 keV. 

Figure [2] shows that the RMS and area pulsed fractions 
measureafrom Chandra HRC between 0.5-10 keV are 
higher than the results measured from NuSTAR in 3-10 
keV, however, the Chandra results are consistent with 
the NuSTAR results above 19 keV to within 0.2cr. This 
suggests that the variation of pulsed fraction from NuS¬ 
TAR is due to the variation in the PWN contribution, 
which decreases faster than the pulsar spectrum at 
higher energy. In Section [4| we will use these results to 
estimate contributions from the PWN and the pulsar’s 
off-pulse component, and discuss constraints on the 
photon index of the pulsar’s off-pulse component. 


^ http://heasarc.gsfc.nasa.gov/docs/software/tools/pimms.html 
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3.2. Phase-averaged, and Pulsed Spectral Analysis 

In order to examine the spectrum of the pulsed emis¬ 
sion from the pulsar, we first defined phase interval 0 - 
0.7 as the pulsed component and the mean count rate 
between phases 0 .7-1.0 as the off-pulse component (see 
FigJTb . We extracted 3-79 keV pulsar events and pro¬ 
duced pulsed spectra. The off-pulse level was then sub¬ 
tracted from the pulsed component to yield the pulsed 
spectrum. In every observation, the two telescopes each 
produced a spectrum, so we had eight spectra from these 
four NuSTAR observations (Section^. Each of the eight 
spectra was grouped into at least 600 counts per energy 
bin and then fitted jointly to various models, using the 
software package XSPEC vl2.8.1g and using y - 2 mini¬ 
mization. To check the results, we also grouped the spec¬ 
tra into 20, 50, 100, 200 and 400 counts per bin. In each 
case, we used both the y 2 minimization, and the cstat 
statistic, which is for Poisson data. The results did not 
differ significantly, so we only report the results for the 
spectra that were grouped into 600 counts per bin, and 
fitted with y 2 minimization statistics. 

We first fitted the pulsed spectra with an absorbed 
single power-law model over 3-79 keV. Since the NuS¬ 
TAR data are not sensitive to a small change in the 
equivalent hydrogen column Nn, we fixed it at the pre¬ 
viously reported value 0.95 x 10 22 cm -2 (Gaensler et 
al. 2002). Cross normalization factors were introduced 
for each spectrum because each observation sampled the 
source in a different detector region, as shown in Table 
[2j Two of the normalization factors were below 1, since 
part of the pulsar field fell into the detectors’ chip gaps. 
We then tied all the fitting parameters except the cross 
normalization fact ors fo r each spectrum. Table |3] shows 
the result. Figure |3(a)| shows the spectra and the fitted 
single power-law model. Since we used phase 0-0.7 to 
produce the pulsed spectra, the exposure time was 0.7 
times smaller than it should be. We corrected the nor¬ 
malization parameter K and flux results throughout the 
spectral analysis with a factor of 1/0.7. With y 2 per de¬ 
gree of freedom (dof) 278/254, the null hypothesis that 
the power-law model describes our data had a probability 
of p = 14.41% so t he result was acceptable. However, 
note in Figure 3(a)| that a systematically curved trend 
is seen in the residuals. Curved trends are also visible 
in the residuals of the spectra grouped into at least 400 
counts per bin, and of the rebinned plots of the spectra 
grouped into at least 200, 100, 50 and 20 counts per bin. 
These imply that the single power-law model may not be 
the optimal representation of our data. 


Following the work of Cusumano et al. (2001), we next 
fitted the spectra with a logarithmic parabolic (logpar) 
model: 

( 1 ) 


A(E) = K'iE/Ep)- 1 *- 13109 ^/^. 


Here, a is the slope at the pivot energy ( E p ), f3 is the 
curvature term, K is the normalization factor and we 


fixed E p at 1 keV, as Cusumano et al. (20011 did. The 


values obtained were a = 1.16±0.05, fi = O.iliO.01 and 
K = (2.81^q Jg) x 10 -3 ph cm -2 s _1 keV -1 . A positive /3 
value means that the spectru m bec omes softer at higher 
energy. Table [3] and Figure 3(b) show the results and 
spectra. With y 2 per dof 254/253, the null hypothesis 
that the logpar model describes our data had a probabil¬ 


ity of p= 47.64%, and the model was acceptable. 

To determine whether the improvement in y 2 using the 
logpar model is statistically valid, we performed the stan¬ 
dard F-test to estimate the probability that the improve¬ 
ment in y 2 tests was by chance. We found the probability 
to be 1.5 x 10 -6 , so the logpar model is statistically a bet¬ 
ter description of the da ta than the power law. Also the 
residuals in Figure 3(b) do not show an obvious trend. To 
test the probability that the residuals were random, we 
perf ormed the Runs test (also known as Wald-Wolfowitz 
test; |Orlandini et al.|2012| |Barlow|1989|) in XSPEC. The 
result was 2 % tor the power-law model, so the curved 


trend seen in the residuals of Figure 3(a) is statistically 
significant. For the logpar model, the r esult was 8 %, so 
statistically the residuals of Figure 3(b) are 4 times more 
likely to be randomly distributed than those of Figure 
3(a) Hence, the logpar model provides a statistically 


superior fit to the spectrum compared with the standard 
power law. 

3.3. Broadband Spectrum 

We also produced a 3 keV- 500 MeV broadband spec¬ 
trum for PSR B1509—58 and show it in Figure [4] The 
X-ray data are from NuSTAR and the y -ray data are 


from observations made with CO MPTEL (Kuiper et al. 
1999| ) and AGILE ( Pilia et al.|2010 1 . Though significant 


pulsations were also detected by Fermi LAT (30- 1000 
MeV) (den Hartog et al. in preparation), we did not use 
the results for our spectral analysis, because the Fermi 
SED of PSR B1509—58 was rather uncertain. Table [4] 
shows results obtained in 3-79 keV, 3 keV-30 MeV and 
3 keV-500 MeV by fitting the logpar model (Equation 
IT} with data from NuSTAR, NuSTAR & COMPTEL, 
and NuSTAR & COMPTEL & AGILE, respectively. All 
fittings were acceptable. Note that we have unfolded 
the spectra in two slightly different methods and show 
the corresponding results in Tables [3] and [4] respectively. 
These two methods are compatible, since in the 3-79 keV 
band they produce results that are consistent within 2 cr 
(the second line in Tables [3] and the first line in Table [T}. 
The first method removes well the effects from the obser¬ 
vations’ Ancillary Response Files (ARFs) and Response 
Matrix Files (RMFs), but required data with response 
files. In the broadband spectral analysis, since the re¬ 
sponse files of the 7 -ray data were unavailable to us, we 
unfolded spectra with a different method assuming that 
the RMFs were perfectly diagonal, which is generally a 
fair approximation. 

As shown in Table [4] spectral parameters obtained in 
these three energy intervals were consistent with each 
other within lcr, hence the energy spectrum from 3 keV 
to 500 MeV can be represented with one smooth logarith¬ 
mic parabolic model. Figure [T] extrapolates the model 
obtained from NuSTAR alone (the red solid line) and its 
lcr ranges (the dashed red lines) to the 7 -ray band. The 
solid blue line is the model obtained from a combined fit 
of all data from 3 keV to 500 MeV (see Table 4) and the 
dashed blue lines are its 1 , 2 an d 3 a ranges. _ 

Table |T| also shows results from Cusumano et al. | (2001 1 
in two energy bands. Results in 2-300 keV were ob- 
tained from BeppoSAX and broadband results in 2 keV 
to 30 MeV were from a combined fit of BeppoSAX and 
COMPTEL. Their broadband spectral parameters a and 
/3 are consistent within 1.5cr with ours (3 keV to 500 
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Table 3 

Pulsed and Phase-averaged Spectral Analysis Results 


Model 


Photon Index 

a 

P 

K 

[10 3 ph cm -2 s -1 keV -1 ] 

Flux (3-79 keV) 
[10 -1 °erg cm -2 s -1 ] 

X 2 /dof 

Null Probability 

Power Law 

1.386 ±0.007 



3.70 ± 0.07 

1.198 ±0.014 

278/254 

14.4% 

Logpar 



1.16 ±0.05 

0.11 ± 0.02 


1.136 ±0.018 

254/253 

47.6% 

Logpai 



0.96 ±0.08 

0.16 ±0.04 

1.8 ± 0.3 

~ 1.33 




a These parameters are from Cusunamo et al. (2001) in the 2- 300 keV band. 


Table 4 

Broadband Spectral Analysis Results at Different Energy Bands 


Energy 

a a 

p a 

K a 

[10 3 ph cm -2 s 1 keV 1 ] 

X 2 /dof 

Null Probability 

NuSTAR Results 

3-79 keV 

1.09 ±0.04 

0.13 ±0.02 

2.47 ± 0.12 

295/260 

7% 

3 keV-30 MeV 

1.104 ± 0.018 

0.124 ± 0.007 

2.51 ± 0.06 

298/263 

7% 

3 keV-500 MeV 

1.083 ±0.015 

0.134 ± 0.005 

2.45 ±0.06 

301/265 

6% 

Cusumano et al. (2001) 

2-300 keV 

0.96 ±0.08 

0.16 ±0.04 

1.8 ±0.3 



2 keV-30 MeV 

1.03 ±0.05 

0.13 ±0.02 

1.9 ±0.1 




a The parameters a, f} and K were defined in Equation^ 

MeV), but their normalization parameter K is ~4.7 cr 
smaller than ours. In this paper, when comparing two 
results, a refers to the root-mean-square of uncertainties 
on the two values unless stated otherwise. 

3.4. Phase-resolved and Pulsed Spectral Analysis 


Ge et al. (2012) obtained phase-resolved X-ray spectra 
ofPSR B1509—58 with the Rossi X-Ray Timing Explorer 
(RXTE) from 3-30 keV and reported variations of the 
phase-resolved photon indices. We used the NuSTAR 
observation to extend this analysis up to 79 keV and 
compare to their results. Moreover, we fitted the spectra 
with both the simple power-law and the logpar models, 
and compared the fits in each phase interval. 

We performed a phase-resolved and pulsed spectral 
analysis for 26 phase bins in the energy range of 3-79 
keV. For this analysis, the pulsar events were extracted 
with a circular region of ape rture radius 60 ,/ centered on 
the pulsa r position ( Sectio n 3.1 ). To compare with the 


results of Ge et al.| 
interval as they die 


we used the same off-pulse 
The* pulse profile in the 3-79 keV 


band is illustrated in the top panel of Figure [5] We de¬ 
fined phase interval 0 to 0.8 as the pulsed portion and 
the mean count rate between phases 0.8 to 1.0 as the off- 
pulse level. The off-pulse level was then subtracted from 
the pulsed portion, yielding the pulsed spectrum. Then 
each pulsed spectrum was divided into 26 phase bins to 
produce phase-resolved and pulsed spectra. Each phase 
interval had at least 5700 counts in 3-79 keV. To compare 
results with Ge et al. (2012; see Fig. [5] middle panel), 
we chose the same phase zero as they aid, and divided 
phase bins using the same phase values. Each spectrum 
was grouped into at least 20 counts per bin and then 
fitted jointly to the power-law model using XSPEC. We 
checked the results with both \ 2 minimization and the 
cstat statistic. The results did not differ significantly, so 
we only report th e res ults from \ 2 statistics. For reasons 
stated in Section 3.2 we fixed iVn at 0.95 x 10 22 cm -2 


and tied all the fitting parameters except the cross nor¬ 
malization factors of the eight spectra in each phase bin. 

The power-law model is acceptable in all phase inter¬ 
vals. Table [5] lists the results and Figure [5] shows the 
photon indices as a function of phase. The photon in¬ 
dices decrease from about 1.53 to 1.34 within the rising 
phase interval, keep stable around 1.34 between phases 
0.28 to 0.44, and start to increase again around phase 
0.54. As shown in the middle and bottom panels of Fig¬ 
ure [5} al l resul ts are consistent within 3tr with those of 
Ge et al. ( |2012| ) obtained by RXTE in the 3-30 keV band. 


i'he logpar model is also acceptable in all phase inter¬ 
vals. Table [5] lists the results and Figure [6] shows the 
parameters as a function of phase. Unlike the photon 
index of the power-law model, the logpar model parame¬ 
ter a does not have an obvious evolving trend within the 
phase interval ~0-0.22, but decreases in the phase in¬ 
terval 0.22-0.3, keeps stable after the phase peak within 
the phase interval of ~0.3-0.44, and shows no obvious 
trend after phase 0.44 (Fig. [6] top panel). The middle 
panel of Figure [6] shows the parameter (3 as a function 
of phase. Note that the (3 values are negative in phase 
intervals 0.22-0.24 and 0.58-0.6, which means that the 
spectra become harder in those two phase intervals. In 
the other 24 phase intervals, /? values are positive, which 
implies that the spectra become softer. A larger positive 
/3 value indicates that the spectra become soft faster in 
that phase interval, but our parameter [3 does not show 
an obvious evolving trend with phase. 

To compare the power-law model to the logpar model 
at different phase intervals, we listed their null hypothe¬ 
sis probabilities obtained from the \ 2 minimization tests 
in Table [5] and plotted them in the bottom panel of Fig¬ 
ured)] Both models were acceptable in all phase intervals, 
andin 20 out of the 26 phase intervals the null hypothesis 
probabilities were above 68%. The logpar model did not 
make a significant improvement compared to the power 
law model in any of the 26 phase intervals. The differ- 
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Figure 3. (a) Top panel: NuSTAR observation data and the 

power-law model. Bottom panel: residual (data — model). 
The power-law model is consistent with our data, y 2 P er dof is 
277.989/254 (p = 14.41%). Note that a systematically curved 
trend is seen in the residuals, (b) Top panel: NuSTAR obser¬ 
vation data and the logpar model. Bottom panel: residual (data — 
model). The logpar model fits with our data better. y 2 per dof is 
253.662/253 (p = 47.64%). No visible trend is seen in the residuals. 
Eight Colors correspond to the eight spectra used in this analysis. 



Figure 4. The pulsed broadband energy spectrum of PSR 
B1509—58 from soft X-ray to 7 -ray. The red squares are the AG¬ 
ILE observa tions in the 30-100 MeV and 100-500 MeV bands (|Pilia| 
let al.| 20T5t. The arrow is the AGILE upper limit above 500 MeV 
band fpilia et al.||2010|). The blue circles are the COMPTEL ob¬ 
servations m trie U. Yb—3 MeV, 3-10 MeV and 10—30 MeV bands 
(|Kuiper et al.||1999|). Data from 3-79 keV are from our NuSTAR 
observations. The solid red line is the logarithmic parabolic model 
obtained from NuSTAR data, and the dashed red lines show its lcr 
range. The solid blue line is the logpar model obtained from all 
data from 3 keV to 500 MeV, and the dashed blue lines are its 1, 2 
and 3cr confidence regions. 



Phase 


ences between the null hypothesis probabilities of the two 
models were less than 4% in 25 out of the 26 phase in¬ 
tervals, and the largest improvement in null hypothesis 
probability was 6.4%, which is insignificant. In addition, 
the bottom panel of Figure [6] illustrates that the null hy¬ 
pothesis probabilities for neither model has any obvious 
correlation with the phase values. 

4. DISCUSSION 

We have reported on a pulsed spectral analysis and 
timing analysis of PSR B1509—58 with observations 
made by NuSTAR and Chandra , as well as previous re¬ 
sults measured by COMPTEL and AGILE. The phase- 
averaged and pulsed spectral analysis shows that a 


Figure 5. Phase-resolved and pulsed spectral fits with the power 
law model. Top panel: photon indices of PSR B1509—58 at dif¬ 
ferent phase ranges (black filled circles). The green line illustrates 
the pulse profile. Middle: comparison of the photon indices ob¬ 
tained in this paper (open black circles) with those of Ge et al (red 
filled diamonds). Bottom: di fference between the photon indices 
obtained in this paper and in |Ge et ah] ( | 2012 | ). 


smooth logarithmic parabolic model represents our 3-79 
keV NuSTAR data better than a simple power-law model 
(Fig. | and Table [3]). In addition, the logpar model de¬ 
scribes our broadband spectrum well from 3 keV up to 
500 MeV (Fig. [ 5 ] and Table [IJ. The phase-resolved and 
pulsed spectral analysis shows that the photon indices in 
26 ph ase int ervals are consistent with the results of |Ge| 
et al. (2012) (Fig. [5] and Table [5j. However, the logpar 
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Table 5 

The Phase-resolved and Pulsed Spectral Analysis Results 


Phase 

Power Law 

Logpar 

Photon Index 

x 2 

dof Null Hypothesis 
Probability 

a 

P 

x 2 

dof Null Hypothesis 
Probability 

0.00-0.14 

1.530±0.123 

1207.66 

1126 

0.045 

1.026T0.852 

0.250 ±0.418 

1207.35 

1125 

0.044 

0.14-0.16 

1.344±0.064 

264.861 

285 

0.798 

0.615 ±0.488 

0.357 ± 0.235 

262.49 

284 

0.815 

0.16-0.18 

1.560±0.046 

297.938 

352 

0.983 

1.273±0.318 

0.145 ±0.158 

297.09 

351 

0.983 

0.18-0.20 

1.456±0.034 

417.727 

436 

0.727 

1.251 ± 0.232 

0.100 ±0.112 

416.925 

435 

0.725 

0.20-0.22 

1.424±0.026 

573.546 

550 

0.236 

1.121 ± 0.182 

0.149 ±0.089 

570.705 

549 

0.252 

0.22-0.24 

1.456±0.022 

634.603 

682 

0.902 

1.515 ±0.136 

-0.029 ±0.066 

634.411 

681 

0.899 

0.24-0.26 

1.445±0.018 

720.945 

805 

0.984 

1.248 ±0.119 

0.096 ± 0.058 

718.06 

804 

0.986 

0.26-0.28 

1.398±0.017 

874.19 

872 

0.473 

1.164 ± 0.109 

0.113 ±0.052 

869.266 

871 

0.510 

0.28-0.30 

1.366±0.016 

813.851 

905 

0.986 

1.050 ±0.108 

0.153 ±0.052 

804.684 

904 

0.992 

0.30-0.32 

1.358±0.016 

843.144 

889 

0.862 

0.948 ±0.109 

0.196±0.052 

827.622 

888 

0.926 

0.32-0.34 

1.367±0.017 

825.724 

874 

0.877 

1.184±0.110 

0.088 ± 0.053 

822.84 

873 

0.886 

0.34-0.36 

1.327±0.018 

741.81 

818 

0.973 

0.951 ± 0.119 

0.181±0.057 

730.939 

817 

0.986 

0.36-0.38 

1.350±0.018 

669.523 

812 

1.000 

1.091 ± 0.117 

0.123 ±0.055 

664.362 

811 

1.000 

0.38-0.40 

1.341±0.018 

728.888 

779 

0.900 

0.983 ±0.121 

0.172 ±0.058 

719.396 

778 

0.934 

0.40-0.42 

1.341±0.019 

735.784 

774 

0.834 

1.106 ±0.125 

0.113 ±0.060 

732.077 

773 

0.851 

0.42-0.44 

1.364±0.019 

747.241 

758 

0.603 

1.153 ±0.130 

0.103 ±0.062 

744.481 

757 

0.620 

0.44-0.46 

1.388±0.020 

718.338 

725 

0.563 

1.316 ±0.132 

0.035 ±0.063 

718.03 

724 

0.556 

0.46-0.48 

1.382±0.021 

635.89 

676 

0.863 

0.988 ±0.145 

0.191 ± 0.070 

628.223 

675 

0.901 

0.48-0.50 

1.400±0.023 

618.462 

617 

0.476 

1.276 ±0.154 

0.060 ± 0.074 

617.793 

616 

0.472 

0.50-0.52 

1.397±0.025 

516.611 

565 

0.928 

0.957 ± 0.176 

0.215 ±0.086 

510.051 

564 

0.949 

0.52-0.54 

1.461±0.028 

459.932 

513 

0.955 

1.126±0.189 

0.163 ±0.092 

456.378 

512 

0.963 

0.54-0.56 

1.410±0.032 

483.215 

464 

0.260 

1.040 ±0.216 

0.180 ±0.104 

480.017 

463 

0.283 

0.56-0.58 

1.432±0.037 

385.549 

410 

0.802 

0.913 ±0.267 

0.258 ±0.131 

381.653 

409 

0.830 

0.58-0.60 

1.442±0.045 

322.192 

361 

0.930 

1.581 ± 0.294 

-0.069 ±0.142 

321.971 

360 

0.926 

0.60-0.64 

1.495±0.038 

565.61 

614 

0.919 

1.055 ±0.262 

0.218 ±0.128 

562.616 

613 

0.928 

0.64-0.80 

1.549±0.056 

1306.21 

1424 

0.988 

0.972 ± 0.388 

0.288 ±0.192 

1303.88 

1423 

0.989 



Phase 

Figure 6. Phase-resolved and pulsed spectral fits with the log- 
par model (Equation 1). The black line in each panels shows the 
pulse profile. Top panel: The parameter a as a function of phase. 
Middle: The parameter f3 as a function of phase. Bottom: Null 
hypothesis probabilities obtained from x 2 tests for the power-law 
model (red open squares) and the logpar model (blue open stars). 
Many of them are very close to each other. 


model does not describe our spectra significantly better 
than the power-law model in any of these phase intervals 
(Fig. [6]). In the timing analysis, we found that the pulsed 
profiles in 14 energy bands between 3-79 keV show no 
significant differences (Fig. III. We measured the pulsed 
fraction of PSR B1509—58 m the hard X-ray band for 
the first time (Fig. [2]). The measured pulsed fraction in¬ 


creased from ~0.68 to ^0.88 from 3 to 19 keV because of 
the PWN contamination which is significant due to the 
broad NuSTAR PSF, and then reached a plateau near 
the latter value within 19-79 keV. 


4.1. Pulsed and Phase-Averaged Spectrum 


In comparison to the results of Cusumano et al. (2001) 
(Table [3j, uncertainties in our parameters are 50%-75% 
smaller. To compare the shape of our logpar model 
to Cusumano’s, we fixed all parameters at Cusumano’s 
values but let the normalization factors vary, and fit¬ 
ted our NuSTAR data with this model in XSPEC with 
X 2 minimization. The results were not acceptable, with 
a probability smaller than 10” 8 . Figure [ 7 ] compares 
Cusumano’s model to our logpar model, and to the NuS¬ 
TAR data. The difference in shape is visible. More¬ 
over, our pulsed flux at 2-10 keV is (2.39 ±0.04) x 10” 11 
erg cm” 2 s , significantl y larger than the 2.0 x IQ” 11 
erg cm” 2 s” 1 reported by Cusumano et al. (20011, but 
consistent w ith the (2 .30 ± 0.02) x 10 _ii erg cm” 2 s” 1 
reported by Ge et al. (20121, and the (2.52 ± 0.05) x 
10~ 11 e rg cm” 2 s” 1 reported by Kuiper fe Hermsen| 
(2015). However, Cusumano’s results were determined 
from BeppoSAX observations. As argued in Section JT| 
we consider our NuSTAR results more reliable than the 
BeppoSAX results, as BeppoSAX had systematic cross¬ 
calibration discrepancies among its three instruments 
among different energy bands. 

We also compared the results of our observations with 
previous X-ray observations and summarize the compar¬ 
ison in Table [6] A simple extrapolation of the observa¬ 
tion results of OSSE (50 keV-5 MeV), Ginga (3-60 keV), 
































































Fig ure 7. Energy spectrum of our NuSTAR observation and that 
from|Cusumano et al.|(|2001|). The g reen points are NuSTAR data. 
The solid black curve is the result of|Cusumano et al.| ([2001 ). The 
solid red curve is the NuSTAR logpar spectral model, and the 
dotted curves are its 1, 2 and 3a ranges. 


Table 6 

Predicted pulsed flux 


Observation 

Energy band 
[keV] 

Photon Index 

Flux in 3-79 keV 
[IQ-iOgrg cm -2 s -1 ] 

NuSTAR 

3-79 

1.386 ±0.007 

1.136 ±0.018 

Ginga 

3-60 

1.30 ±0.05 

~1.50 

RXTE 

2-250 

1.358 ±0.014 

~2.20 

WELCOME 

94-240 

1 64 +u ' 4;:5 

± -°^-0.42 

~3.14 

OSSE 

50-5000 

1.68 ±0.09 

~1.26 


RXTE (2-250 keV) and WELCOME (94-240 keV) into 
our energy range predicts a pulsed flux of ~1.26xlO -10 , 
~1.50xl0“ 10 , ~2.20xKT 10 and ~3.14xlO” 10 erg cm" 2 
s -1 , respectively. The pulsed fluxes predicted by OSSE 
and Ginga are consistent with our NuSTAR flux in 3-79 
keV, while RXTE and WELCOME each predict a pulsed 
flux of ~2 and ^3 times larger than our NuSTAR flux. 
Note that cross-calibration discrepancies among NuS¬ 
TAR, Ginga and RXTE against the Crab nebula have 
been found. NuSTAR was reported to measure the Crab 


spectrum ~15% lower than G inga and RXTE (Madsen 
et al.|2Qj~5 Kirsch et al.|2005 ). Also NuSTAR has an un¬ 


certainty of ~5% on its absolute flux measurement, and 
the C rab fl ux could chang e by a few percent per year 
(Wilson-Hodge et al.|2011). Furthermore, the Ginga (3- 
6(TkeV')7T?XT , E r ^"2-250 keV), WELCOME (94-240 keV) 
and OSSE (50 keV-5 MeV) observations were reported to 
fit with simple power-la w models. Their photon indices 
were 1.3 0T0.05 ( Ginga; |Kawai et al.|1993| ), 1.358T0.014 
(RXTE; |Marsden et al.||1997|), 1.64+^ ( WELCOME; 


Gunji et al. 1994|) and 1.68 ± 0.09 (OSSE; Matz et al. 


1994p (Table [b|). This indicates that the pulsed spec¬ 


trum becomes softer as energy increases, consistent with 
the behavior of our logpar model. 

Figure [8] shows a comparison of the best-fit energy 
spectrum from our NuSTAR observations with those of 
the four previous X-ray observations. OSSE is within 3a 



Figure 8. Energy spectrum of NuSTAR data and previous X-ray 
observations. The solid red line is the NuSTAR logpar spectral 
model, and the dotted red lines are its 1, 2 and 3 a ranges. The 
solid black line, solid cyan line, solid blue line, and solid magenta 
line, from top to bottom, are X-ray spectra observed by RXTE , 
WELCOME , Ginga and OSSE , respectively. 

from our NuSTAR spectrum, while RXTE , WELCOME 
and Ginga do not agree well with our new NuSTAR mea¬ 
surements. This perhaps is due to the cross-calibration 
issues between instruments described earlier in this sec¬ 
tion, since t he source flu x has be en ve ry stable over a 
long period (Livingstone & Kaspi 2011). Note that the 
uncertainty of the model increases dramatically with en¬ 
ergy. 

In t heir broadband spectral analysis, |Cusumano et ah] 
(20011 claimed that a smooth logarithmic parabolic 
model could represent the spectrum well up to 30 MeV, 
but must break more rapidly afterwards. However, our 
result shows that a logarithmic parabolic function can 
represent the spectrum well up to 500 MeV. Moreover, 
one interesting feature of the broadband energy spectrum 
is its peak energy. According to Equation H] the energy 
flux should have a maximum at the following energy: 


E m = E 0 10 [(2 -“ )/(2/3)] , 


( 2 ) 


which corresponds to 3 MeV and a 90% upper limit of 8 
MeV, using the best-fit parameters obtained from NuS¬ 
TAR alone, 4±2 MeV using the combined fit of NuSTAR. 
and COMPTEL from 3 keV to 30 MeV, and 2.6T0.8 MeV 
using results obtained from all data up to 500 MeV. This 
value is consist ent with the recent results of Kuiper & 
Hermsen (20151. Using multi-year data across ~2 keV-1 
GeV from RXTE PCA/HEXTE , Fermi LAT , COMP¬ 
TEL and INTEGRAL ISGRI, they performed a broad¬ 
band spectral fit of the pulsed emission of PSR B1509—58 
to a power-law model with a modified exponential cutoff, 
and derived a maximum luminosity p er en ergy d ecade at 
^2.5 MeV. The result of Cusumano et al. (2001|) was ob¬ 
tained from a combined fit of their BeppoSAX X-ray ob¬ 
servation and COMPTEL, from 1 keV to 30 MeV. Their 
peak energy was 5 MeV, with a 90% upper limit of 14 
MeV. Their value agrees with all of ours, though their un¬ 
certainty is ^9 times larger than what we have obtained 
from the combined fit. Note that each of the y-ray data 
points plays an important role in our combined fits, and 
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in the determination of the peak value. 


Cusumano et al. (20011 also discussed the relation of 
the peak values between PSR B1509—58 and the Crab 
pulsar. For the outer gap model, the synchrotron spec¬ 
trum emitted from the secondary electron-positron pairs 
created by inward flowing curvature photons of Crab-like 
pulsars is predicted to o bey the following scaling relation 
(Cusumano et al. 200 l| ): 

( 3 ) 


(■ E m /E mC ) ~ ( P/P c ) 7 l \P/Pc ) 1/4 - 


Here, Pc and Pc are period and spin-down rate of the 
Crab, E m , E m c are energies corresponding to the spec¬ 
tral peaks of the Crab-like pulsar and the Crab. It is 
suggested that a Crab-like pulsar spectrum that follows 
this relation is likely to be synchrotron rad iation from 
the secondary pairs ( Cusumano et al. 20011. For PSR 
B1509—58, the RHS of Equation [3] is ~20. For the 
Crab, the SED peak energy ( E m c ) is ~14 keV for the 
main pulse, and ~200 k eV for the interpulse component. 


Cusumano et al. (2001) concluded that their observation 
agreed with this ratio when the spectrum of the Crab’s 
interpulse region was considered, but it was about one 
order of magnitude larger than the result when the spec¬ 
trum of the Crab’s main pulse was considered. We can 
improve on this result thanks to our more precisely mea¬ 
sured peak energy. From our analysis, the peak value 
obtained from our combined fit (3 keV-500 MeV) gives 
a ratio of ~5 ± 2 when considering the Crab’s interpulse 
region, and ~80 ± 20 for the Crab’s main pulse. Neither 
agrees well with the ratio of ~20 predicted by Equation 
[3j Note that the exact uncertainties of E m c or the RHS 
of Equation [3] are unavailable to us, so the conclusion is 
based on estimations. 

4.2. Predicted Spectra from Theoretical Models 

PSR B1509—58 has an atypical spectrum in the X-ray 
and 7 -ray bands. 7 -ray pulsars typically have s pectral 
peaks in the GeV energy band (Thompson 20041, while 
PSR B1509—58 was observed to emit very little 7 -ray 
emission, and its energy spectrum pe aks in the hard X- 
ray ba nd at 2.6 ± 0.8 MeV (Section |3.3|). |Wang et al. 
(|2013I proposed a model in the context of outer gap ac¬ 
celeration to explain the weak 7 -ray emission from this 
pulsar and calculated the X-ray spectrum. They pro¬ 
posed that the outgoing GeV-band 7 -ray radiated via 
curvature radiation from the electron-positron pairs cre¬ 
ated in the null charge surface are missed because of the 
Earth viewing angle, while the inward propagating 7 - 
rays passing near the pulsar surface may be absorbed by 
the magnetic field and converted into electron-positron 
pairs. The secondary pairs that travel with sufficiently 
large pitch angles with respect to the magnetic field lines 
emit synchrotron radiation which peaks at ~10 MeV, and 
covers a wider sky area compared with curvature radia¬ 
tion of the primary particles. They concluded that only 
synchrotron radiation could be observed due to our view¬ 
ing angles. However, the model may undenjaredict emis¬ 
sions in the X-ray band. For example, Wang et al. (2013) 
compared the energy spectrum predicted by their model 
with previous X-ray observations by Ginga , OSSE , WEL¬ 
COME and RXTE, and found that the model predicts 
lower energy flux than is seen in all these observations 
except that of OSSE. The predicted energy flux is ~half 


as much as observed by RXTE and WELCOME , and 
more than 5<r lower than the RXTE observation. 

We quantitatively compared our results with the the- 
or etical model prediction that was shown in Figure 6 of 
Wang et al. (2013). By inspection, their theoretical pre¬ 
diction is roughly consistent with our NuSTAR 3ct ranges 
(Fig. [ 8 |. Mo reover, the theore tical model predicts a peak 
at ~3 MeV ( Wang et al.|2013 1, which agrees with our re¬ 
sults. 

4.3. Pulsed and Phase-resolved Spectrum 

We also studied the pulsed, phase-resolved spectrum 
of PSR B1509—58. We found that the photon index is 
1.34±0.02 in the middle of the pulse, becoming larger at 
the two sides, for the power-law model. O ur resu lts are 
largely consistent with those of Ge et al.] 12012[ ). Note 
that in 24 ph ases our photon indices are sli 
than those of Ge et al. ( 2012| ) (see Fig. [5 


lilghfiy larger 
middle and 

bottom panels'^ This may be consistenCwith the fact 
that the photon index increases with energy in the X- 
ray band, since our results were obtained from NuSTAR. 
between 3-79 keV, higher than their 3-30 keV energy 
band. 


Ge et al. (2012) performed pulsed and phase-resolved 
spectral analyses for three young and bright X-ray 
pulsars: the Crab pulsar, PSR B1509—58 and PSR 
B0540—69. Their results show that the photon indices 
of PSRs B1509—58 and B0540—69 are harder at the cen¬ 
tres of the pulse profiles and softer at the wings, while 
those of the Crab pulsar are softer at its two peaks. How¬ 
ever, their timing analysis shows that the pulse profiles of 
PSRs B1509—58 and B0540—69 could each be further de¬ 
composed into two narrower Gaussian components. Tak¬ 
ing this into account, they concluded that overall the 
three pulsars’ photon indices had a similar evolving trend 
with their pulsed flux. Our NuSTAR observations sup¬ 
port the same conclusion. 

4.4. Pulsed Fraction 


Wang et al. 2013 )’s model clearly predicted X-ray off- 
pulse emission from the pulsar, and those authors sug¬ 
gested that it was probably produced by electrons accel¬ 
erated transverse to the field line. We have measured the 
pulsed fraction of PSR B1509—58 with NuSTAR in the 
3-79 keV band and with Chandra HRC in the 0.5-10 keV 
band. The observed off-pulse emission originated from 
both the pulsar and the PWN. With our pulsed fraction 
measurements, we can determine the contributions from 
the PWN and the PWN-free pulsar off-pulse components 
to the total emission, and estimate the photon index of 
the pulsar’s off-pulse component for the first time. 

Based on the results shown in Figure [2j we had that 
at 3-4 keV, the PWN-free pulsed fraction measured with 
Chandra is: 


PP 


PWN-free — 


pulsed 


pulsed + o//pidsepwN— free 


= 0.88±0.06, 


while the PWN-included result from NuSTAR is: 

pulsed 


PFw itliPWN = 


pulsed + offpulse PWN _ free + PWN 
= 0.68 ±0.04. 

From these two results, we derived that the PWN 
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contrioff-pulsebutes 0.23 ± 0.06 of the total number of 
NuSTAR photons, and that the pulsar’s PWN-free off- 
pulse component contributes 0.09±0.07 of the total num¬ 
ber of photons. The latter result is roughly consistent 
with what is shown in Figure [2] where the pulsed frac¬ 
tion approaches 90%. The number of photons detected 
from the PWN is 3 ± 2 times of those from the pul¬ 
sar’s PWN-free off-pulse component at ~3-4 keV. Con¬ 
sequently, the pulse-off components at soft X-ray bands 
are likely to be PWN-dominated, and as energy increases 
the PWN contribution decreases until ~19 keV, where 
it no longer matters. The off-pulse components of 19- 
79 keV are composed of the pulsar’s PWN-free off-pulse 
component alone. 

Additionally, we can estimate the photon index of the 
pulsar’s PWN-free off-pulse component. The PWN spec¬ 
trum is much softer than the pulsar’s: th e photo n i ndices 
of the PWN and the pulsar are 1.6-1.8 (An et al.|20141 
and 1.386 ± 0.007, respectively. The P WIN is also much 
fainter than the pulsar: at 19 keV, its flux is less than 
1% of the pulsar’s. Hence we can safely neglect the PWN 
contribution above 19 keV and regard the NuSTAR re¬ 
sults between 19-25 keV and 25-79 keV as PWN-free 
pulsed fractions of the pulsar. Meanwhile, the Chan¬ 
dra HRC results also give the PWN-free pulsed fraction. 
Figure [2] shows that the four PWN-free pulsed fractions, 
namely the two NuSTAR results in 19-25 keV and 25- 
79 keV, as well as the two Chandra results within 0.5-10 
keV, are consistent. Assuming that the PWN-free off- 
pulse components is well represented by a simple power- 
law model, the PWN-free pulsed fraction is then 


PF = 


I\E 


-a-/3Zog(#-) 


KE-a-piog^) + Koff _ pulseE - r 


( 4 ) 


Here K 0 f / -pulse and F are the normalization factor and 
the photon index of the pulsar’s PWN-free off-pulse com¬ 
ponent, to be determined. AT, a and /? are the same as 
in Equation and Table [3j With these four PWN-free 
pulsed fraction results, we took their upper and lower 
la limits and estimated that the off-pulse component’s 
photon index is between 1.26 and 1.96. This constraint 
is weak due to the large uncertainties. NuSTAR obser¬ 
vations with longer exposure time could be useful for 
constraining the off-pulse spectrum better. 


5. SUMMARY 

We have presented NuSTAR hard X-ray observa¬ 
tions of the young rotation-powered radio pulsar PSR 
B1509—58 in the supernova remnant MSH 15— 52. We 
have confi rmed the curvature in t he hard X-ray band re¬ 
ported by Cusumano et al. ([2001), and have showed that 
the logpar model is statistically a better representation 
of our spectrum in both the NuSTAR. energy band from 
3 to 79 keV, and a broadband from 3 keV to 500 MeV. 
The logpar model predicts that the broadband energy 
spectrum peaks at 2.6 ± 0.8 MeV, ~9 times more pre¬ 
cisely determined than befo re, and cons istent with the 
recent result of Kuiper & Hermsen (2015). We have also 
fitted the logpar model and the power-law model to the 
phase-resolved spectra in 26 phase intervals, but found 
that statistically the logpar model was not significantly 
superior to the power-law model in any of these phase 


intervals. In addition, we have measured the pulsed frac¬ 
tion of PSR B1509—58 in the hard X-ray energy band 
for the first time. The results imply that the PWN con¬ 
tributes 0.23 ±0.06 of the total number of NuSTAR. pho¬ 
tons, while the pulsar’s PWN-free off-pulse component 
contributes 0.09 ± 0.07 of the total number of photons. 
We have also estimated that the off-pulse component’s 
photon index is between 1.26 and 1.96. To improve the 
constraint, NuSTAR observations with longer exposure 
time could be helpful. Fina lly, these new NuS TAR. ob¬ 
servations support a model (|Wang et al.||2013| in which 
the X-rays originate from secondary pairs producing syn¬ 
chrotron radiation. The model describes how outgoing 
GeV emission goes undetected because of the Earth view¬ 
ing angle. This helps to understand PSR B1509—58’s 
lack of GeV emission, and provides a possible answer 
to the long-lasting question about the source’s atypical 
spectrum. 
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